Ab initio study of charge order in Fe304 
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We present a self-interaction corrected local spin density (SIC-LSD) study of the electronic struc- 
ture and possible charge order of magnetite, Fe304. The issue of charge order in magnetite is 
explored in both cubic and orthorhombic structures, the latter being an approximation to the true, 
low temperature, monoclinic structure. We find that the Verwey charge ordered phase is not the 
groundstate solution for this compound either in cubic or orthorhombic structure. We conclude that 
the structural distortions, more than localization/delocalization correlations, are responsible for the 
charge disproportionation in the low temperature phase. 



I. INTRODUCTION 

Based upon its high magnetoresistive properties, mag- 
netite, Fe304, is a system of interest for such technolog- 
ical apphcations as computer memory, magnetic record- 
ing, and other "spintronic" devices. It has been exten- 
sively studied for the past 60 yearai, but despite such an 
extended effort the issue of a possible charge order below 
the so-called Verwey transition temperature, Tv=122 K, 
has not been understood yeli^. Above the Verwey tran- 
sition temperature, magnetite is thought to be a half- 
metal with a highest known Tj, of 860 K. At the tran- 
sition temperature, the conductivity of the system de- 
creases abruptly by two orders of magnitude, as a small 
gap opens up in the density of states spectrum.^ This 
first order metal-insulator transition'* is associated with 
a lattice distortion from cubic to monoclinic, which how- 
ever has not been fully resolved.-"-'- The transition has 
been viewed as an order-disorder transition in relation to 
the arrangement of cations on the octahedral sites of the 
inverse spinel structure whose formal chemical formula 
can be written as Fe^'^"'' [ Fe^+, Fe^+ ]b04^~, where 
A and B refer to the interstitial, i.e., tetrahedral and 
octahedral sites, respectively. This simple ionic picture 
implies that magnetite is a mixed-valent compound. The 
A sublattice is occupied only by Fe'^+ cations in a ferro- 
magnetic arrangement, while the B sublattice is occupied 
by both Fe^"*" (Bl-site) and Fe^+ (B2-site) cations, also 
in a ferromagnetic manner which, however, is antiferro- 
magnetically aligned to the ferromagnetic order on the A 
sublattice. The oxygen atoms occupy the face-centered- 
cubic (fee) close-packed sites. The Fe^+ cation can be 
viewed as Fe'^"'' plus an "extra" electron, hopping freely 
between the octahedral sites above the Verwey transi- 
tion, while below the transition these "extra" electrons 
become localized or quasi-localized and some kind of long 
range order sets in on the B sublattice. Verwey, who 
argued that below the transition temperature, Ty, the 
Fe"^"*" and Fe^"*" cations order in alternate (001) planes, 
looked upon this transition as an electron localization- 
delocalization transition. The Mossbauer spectra clearly 
demonstrate existence of two types of localized ions on 



the octahedral sitesn^ The low temperature structure is 
a rhombohedral distortion of the cubic spinel arrange- 
ment to the first approximatiouf^i but this splits in the 
ratio 3:1 for the occurance of the Fe^"*" and Fe'^"'" cations 
on the corner-sharing tetrahedra that are formed by the 
octahedral sites. This is not consistent with the Ander- 
son condition of minimum electrostatic repulsion, imply- 
ing the ratio of 2:2.12, A recent combined high resolution 
X-ray and neutron powder diffraction study*" confirmed 
a small charge disproportionation of about 0.2 electron 
on the octahedral sites, and a violation of the Anderson 
condition. The current experimental situation suggests 
that the long range order present below the Verwey tran- 
sition is much more complicated than the charge order 
first proposed by Verwey. Some experiments, claim that 
there exists no charge order on the octahedral sites>ii*i^ 

In this paper we present an application of the self- 
interaction corrected local spin density (SIC-LSD) ap- 
proximation to Fe304. The SIC-LSD method allows the 
ab-initio realization of valencies through total energy cal- 
culations for configurations of Fe^+ and Fe^+ cations. 
The method includes all bandstructure effects of the de- 
localized Fe and O electrons, including their hybridiza- 
tion, as well as the electronic structure of the localized 
Fe electrons. As such this methodology has been able 
to describe the valencies of the rare earths successfullyii^ 
Here, we concentrate on the total energies for different 
d electron configurations on the octahedral sites, defin- 
ing Fe valencies, and try to establish whether the simple 
charge order postulated by Verwey is the ground state 
solution. Apart from the Verwey order we also consider 
other obvious arrangements of Fe valencies on the octahe- 
dral sites. There have been other first-principles studies 
of magnetiter^*i^ii& however, this is the first time that 
different Fe valencies and charge order scenarios have 
been studied by an ab initio method. In Refs. andU^ 
the LDA+U method has been used to study the Verwey 
insulating phase albeit in the cubic structure only, and 
neither the total energies nor other charge order phases 
have been considered. In Ref. [Til in addition, the optical 
properties of magnetite have been studied. 

The remainder of the paper is organized as follows. 
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In the next section, we briefly present the theoretical 
background of the SIC-LSD approach. In section III, the 
technical and computational details regarding the appli- 
cation to magnetite are elaborated upon. In section IV, 
we present total energy and density of states calculations. 
We also study the charge order and different Fe valencies 
on the octahedral sites, and try to establish the ground 
state configuration. The paper is concluded in section V. 



II. THEORY 

The basis of the SIC-LSD formalism is a self- 
interaction free total energy functional, E^^'-' , obtained 
by subtracting from the LSD total energy functional, 
E^^^ , a spurious self- interaction of each occupied elec- 
tron state f/'oF^j namely 

occ. 

j^SIC ^j^LSD ^Y^gSIC^ 
a 

Here a numbers the occupied states and the self- 
interaction correction for the state a is 

S'J^ = U[n^]+E^,'^[nc,], (2) 

with [/[ric] being the Hartree energy and E^^^[na] the 
LSD exchange-correlation energy for the corresponding 
charge density Ha and spin density 

The SIC-LSD approach can be viewed as an extension 
of LSD in the sense that the self-interaction correction is 
only finite for spatially localized states, while for Bloch- 
like single-particle states E^^'^ is equal to E^^^ . Thus, 
the LSD minimum is also a local minimum of E^^^ . A 
question now arises, whether there exist other competi- 
tive minima, corresponding to a finite number of localized 
states, which could benefit from the self-interaction term 
without loosing too much of the energy associated with 
band formation. This is often the case for rather well lo- 
calized electrons like the 3 d electrons in transition metal 
oxides or the 4/ electrons in rare earth compounds. It 
follows from minimization of Eq. that within the SIC- 
LSD approach such localized electrons move in a different 
potential than the delocalized valence electrons which re- 
spond to the effective LSD potential. For example, in the 
case of magnetite, five (Fe'^+) or six (Fe^"*") Fe d electrons 
move in the SIC potential, while all other electrons feel 
only the effective LSD potential. Thus, by including an 
explicit energy contribution for an electron to localize, 
the ab-initio SIC-LSD describes both localized and delo- 
calized electrons on an equal footing, leading to a greatly 
improved description of static Coulomb correlation ef- 
fects over the LSD approximation. 

In order to make the connection between valency and 
localization more explicit it is useful to define the nominal 
valency as 

Nval — Z - Ncore - NsiC, 



where Z is the atomic number (26 for Fe), Ncore is the 
number of core (and semi-core) electrons (18 for Fe), and 
Nsic is the number of localized, i.e., self-interaction cor- 
rected, states (either five or six, respectively for Fe'^+ and 
Fe^"*"). Thus, in this formulation the valency is equal to 
the integer number of electrons available for band for- 
mation. To find the valency we assume various atomic 
configurations, consisting of different numbers of local- 
ized states, and minimize the SIC-LSD energy functional 
of Eq. ^ with respect to the number of localized elec- 
trons. The SIC-LSD formalism is governed by the ener- 
getics due to the fact that for each orbital the SIC differ- 
entiates between the energy gain due to hybridization of 
the orbital with the valence bands and the energy gain 
upon localization of the orbital. Whichever wins deter- 
mines if the orbital is part of the valence band or not and 
in this manner also leads to the evaluation of the valency 
of elements involved. The SIC depends on the choice of 
orbitals and its value can differ substantially as a result 
of this. Therefore, one has to be guided by the energet- 
ics in defining the most optimally localized orbitals to 
determine the absolute energy minimum of the SIC-LSD 
energy functional. The advantage of the SIC-LSD for- 
malism is that for such systems as transition metal oxides 
or rare earth compounds the lowest energy solution will 
describe the situation where some single-electron states 
may not be Bloch-like. For magnetite, these would be 
the Fe 3c? states, but not the O 2p states, as trying to 
localize the latter is energetically unfavourable. 

In the present work the SIC-LSD approach, 
has been implemented^^ within the linear muffin- 
tin-orbital (LMTO) atomic sphere approximation 
(ASA) band structure method, in the tight-binding 
representation)^ The electron wave functions are ex- 
panded in terms of the screened muffin-tin orbitals, and 
the minimization of E^^"^ becomes non-linear in the 
expansion coefficients. 



III. CALCULATIONAL DETAILS 

As mentioned earlier, magnetite crystallizes in the 
fee inverse spinel structure with two formula units (14 
atoms) in the primitive unit cell, and the space group is 
Fd3m. In the actual calculations, an additional 16 empty 
spheres have been used in order to obtain better space 
filling and to minimize the ASA overlap. Concerning the 
linear muffin-tin basis functions, we have used 4s, ip, 
and 3d partial waves on all Fe atoms, and treated them 
as low-waves.2^ Including also 4/— basis functions on the 
Fe-sites, treated as intermediate waves, is of no substan- 
tial importance for the final results. On the oxygens only 
2s and 2p partial waves have been treated as low-waves, 
while the 3d— waves have been treated as intermediate. 
On the empty spheres the Is— waves have been consid- 
ered as low- and the 2p's as intermediate-waves. The 
empty spheres are necessary to obtain half-metallic solu- 
tion with the LSD, whose self-consistent charge densities 
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have been used to start SIC-LSD calculations. 

Although, below the Verwey transition, the cubic 
structure undergoes a monoclinic distortion, we have first 
explored the charge ordered phases on the basis of the 
cubic structure and then extended the study to the dis- 
torted structure. Since the true monoclinic structure has 
not been fully resolved, our calculations for the low tem- 
perature phase have been performed for two orthorhom- 
bic structures, given in Refs. iandQ, as closest approx- 
imations to the monoclinic structure. In the latter work, 
the refined atom positions were obtained using a com- 
bination of high resolution X-ray and neutron powder 
diffraction. The respective space groups used are Pmc2i 
and P2/c, each with 56 atoms in the unit cell. How- 
ever, for reliable calculations empty spheres have been 
included, increasing the total number of sites to 104 in 
the first case, and 120 for the refined structure of Wright 
et aU 

To make the SIC-LSD calculations feasible, for the or- 
thorhombic set-ups we had to restrict to minimum basis 
sets, in particular 4s, 4p, Sd— low waves have been used 
on every Fe-atom, 2s,2p— low waves and intermediate 
3d— waves on all oxygen atoms and only Is— low waves 
and intermediate 2p— waves on the empty spheres. Also, 
the number of k-points in the irreducible Brillouin zone 
(BZ) had to be reduced to 18 k-points. 

To study charge order, we have implemented three ob- 
vious scenarios. First, we have considered the simple 
Verwey charge order, where it is assumed that the tetra- 
hedral sites are occupied by Fe'^+ ions, while on all the 
octahedral sites Fe alternates between Fe'^"'" (five d elec- 
trons move in the SIC-LSD potential) and Fe^+ (six d 
electrons move in the SIC-LSD potential) in the suc- 
cessive (001) planes. Subsequently, we refer to this as 
Verwey charge order scenario. In the second scenario, 
we assume that all tetrahedral and octahedral sites are 
trivalent, namely Fe"^^. Finally, in the third scenario we 
assume all the tetrahedral sites to be Fe'^"'", and all the 
octahedral sites to be Fe^+. In what follows, we shall 
refer to these calculations as scenarios 1-3, respectively. 

Note that for a divalent Fe atom one needs to consider 
six d electrons moving in the SIC potential. Five of the 
electrons fill up the full manifold of one spin component, 
whilst the sixth electron has the freedom to populate ei- 
ther of the five states of the other spin component (three 
t2g and two e^ states in the cubic symmetry). As a re- 
sult five separate, self-consistent total energy calculations 
have to be performed, for all different structures and sce- 
narios involving divalent Fe's, to find the lowest energy 
solution. In what follows, only the results of the lowest 
energy solutions are presented and discussed. 



IV. RESULTS AND DISCUSSION 
A. Densities of states 

In Fig. 1 we compare for the cubic structure the LSD 
total and species-decomposed densities of states (DOS) 
with the SIC-LSD counterparts where the Verwey charge 
order has been implemented. In this figure both spin 
components are displayed. In Fig. 2, we compare sep- 
arately the LSD spin-decomposed DOS for the tetrahe- 
dral and octahedral sites with the respective DOS from 
the SIC-LSD calculations, where the octahedral sites are 
additionally split into Fe^"*" and Fe^+. The first thing 
to note about Fig. 1 is that the LSD calculations find 
half-metallic DOS, while a small band gap opens up in 
both spin-channels of the SIC-LSD densities of states, as 
a result of the Verwey charge order implemented on the 
octahedral sites. The gap is about 0.35 eV, which is more 
than two times the experimental gap found in photoemis- 
sion for magnetite below the Verwey transition tempera- 
ture, but in good agreement with the LDA-I-U calculation 
of Anisimov et al.ji^ although the SIC-LSD gap has the 
p-d charge transfer character, while the LDA-I-U gap oc- 
curs between d states of the octahedral Fe ions. Note 
that due to opening of the SIC-LSD gap, the unoccupied 
Fe 3d states have been shifted up in energy, away from 
the Fermi level, while the occupied Fe 3d states appear 
well below the valence band, which is composed mostly of 
oxygen 2p with a small admixture of Fe d states. In the 
LSD plots, the spin-up Fe 3d states occur at the Fermi 
level, while there is a well defined gap in the other spin- 
component. In addition, the predominantly oxygen-like 
valence band is much more structured in the case of LSD, 
due to strong hybridization with Fe 3d states, which ap- 
pear at the top of the valence band. Looking at all differ- 
ent types of Fe sites in both the calculations (see Fig. 2), 
it is clear that the Fe-A sites in the SIC-LSD calculations 
do not contribute to DOS at the Fermi level, so there is 
a well defined gap in both spin channels for these sites. 
The A-sites in the LSD calculation exhibit a very small 
DOS in the spin-up component, and a gap in the other 
spin channel. The octahedral, B-sites, in the LSD calcu- 
lations exhibit a gap in one spin component and in the 
other spin-component well hybridized states can be ob- 
served at the Fermi level, giving rise to the half-metallic 
ground state. Among the octahedral sites in the SIC- 
LSD calculations, the B-sites occupied by Fe'^'*' ions show 
similar DOS to the LSD calculations for the octahedral 
ions, except for the gap appearing in the SIC-LSD case 
in both spin-channels, and the localized spin-down Fe'^"'" 
states occuring well below the valence bands. Regard- 
ing the unoccupied Fe 3d states, the spin-up states for 
the Fe^^ sites lie higher in energy than the correspond- 
ing states for the Fe'^+ sites. The small peak observed 
just below the valence band is due to the sixth localized 
3d state, while the remaining five localized d— states are 
visible at about 8.0 eV in the other spin component. 

In Fig. 3 we show the total and partial Fe- and 
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0-densities of states for the second SIC-LSD scenario, 
namely for the case where all Fe-sites, tetrahedral and oc- 
tahedral alike, are 3+. The DOS are half-metallic, with 
a gap in one spin component and hybridized Fc-0 states 
occuring at the Fermi level in the other spin component. 
The DOS of localized Fe 3d states appear well below the 
valence band at about 13.0 eV. Note substantial differ- 
ences of these DOS with respect to the LSD result of Fig. 
1 . The reason being, that in the LSD all electrons move in 
the effective LSD potential, while in the 3-1- scenario, the 
Fe 3d electrons respond to the SIC potential. Due to this 
localizing potential the total weight of the occupied 3d 
states is shifted well below the valence band, while in the 
LSD picture they appear at the top of the valence band, 
strongly hybridizing with the O 2p states. From the A 
and B-Fe site decomposed spin-polarized DOS in Fig. 4 
we can see that the tetrahedral sites do not contribute at 
the Fermi level. The half-metallic behaviour is entirely 
due to the hybridization of the unoccupied 3d states of 
the octahedral Fe sites with the oxygen 2p states. 

Finally, the LSD and SIC-LSD densities of states for 
the orthorhombic structures (not displayed in this work) 
show similar features to their counterparts in the cubic 
structure, when the respective tetrahedral and octahe- 
dral ions are compared. Specifically, a small energy gap 
of about O.leV and insulating solutions are observed in 
the charge ordered state for both orthorhombic symme- 
tries. Note that in both orthorhombic structures that 
were studied there are six crystallographically different 
B-sites which are subsequently referred to as Bl, B6.^ 
Nevertheless, we still have only two types of Fe ions, diva- 
lent and trivalent. Concerning the tetrahedral (A) sites, 
there are four different crystallographic positions in the 
structure with Pmc2i symmetry, but only two A- type 
positions in the structure with P2/c symmetry. 

Detailed comparison of the densities of states for all 
different crystallographic structures that we have imple- 
mented here, reveals that the second scenario in the or- 
thorhombic structure with the Pmc2i symmetry leads 
to metallic rather than half-metallic DOS. This result is 
caused by the fact that the Fermi level moves slightly 
into the Fe-A 3d peak, which is unoccupied in the cubic 
structure. However, for the refined orthorhombic struc- 
ture of Ref.'', with P2/c symmetry, we again recover a 
half-metallic solution for the all trivalent configurations. 
As in the cubic structure, the Fe-A 3d peak is fully un- 
occupied. This appears to indicate a great sensitivity of 
the half-metallic solution to atom positions in the low 
temperature structure. 



B. Total energies and charge order 

In this subsection we concentrate on the total energies, 
spin magnetic moments, and charge disproportionation 
on octahedral sites, obtained for the scenarios studied 
here. Tables I and II summarize for the cubic and or- 
thorhombic structures, respectively, the spin magnetic 



moments on Fe-sites and energy differences with respect 
to the ground state configuration of the most favourable 
scenario. 

Comparing total energies for the LSD and three dif- 
ferent SIC-LSD charge order scenarios in either cubic 
(Table I) or orthorhombic (Table II) structure, one can 
see that the second scenario, with all Fe's in the triva- 
lent configuration, is the lowest energy solution. It is 
followed by the Verwey charge order scenario, then the 
third scenario, and finally the LSD, which is the most 
unfavourable result. The significant energy difference of 
113 mRy per formula unit (Table I), found between the 
scenario SIC(l) and scenario SIC(2), renders the sim- 
ple Verwey ordering highly unlikely as an explanation 
for the properties in the low temperature phase of mag- 
netite. The last two columns of Table I, AEbf and 
AEs/c representing respectively the loss in band for- 
mation/hybridization energy and gain in localization en- 
ergy, relative to the respective values of the ground state 
configuration, provide further physical understanding of 
this finding. Whilst the localization energy favours more 
localized states, and hence the formation of Fe^+ on the 
B2 octahedral sites, the loss of band formation energy by 
far outweighs the gain in localization energy. So, there 
is a strong competition between localization and band 
formation/hybridization energies, and in the simple Ver- 
wey scenario the loss in band formation energy, due to 
localization of one additional electron on the B2 octahe- 
dral sites, is three times larger than the gain in localiza- 
tion energy. This balance could possibly be reversed by 
increasing the nearest neighbour Fe (B2-type)-0 bond 
length, which could perhaps be realized in a very dis- 
torted structure, although it is not the case for the or- 
thorhombic structures studied here. The expectation is 
that substantially increased Fe (B2-type)-0 bond length 
would lead to an increase in the localization energy of 
all six d electrons of B2-Fe (with 2+ valency), and a re- 
duction of their hybridization, thus establishing the Ver- 
wey phase as the lowest energy solution. However, one 
should keep in mind that we refer to nominal valencies, 
and therefore find that already Fe'^+ contains more than 
five d-electrons (Table HI). All in all, Fe304 turns out to 
have mostly intermediate valencies. 

Inspecting in detail the other numbers of Table I, where 
the results for cubic structure are summarized, one can 
see that because these systems are either insulating or 
half-metallic their total magnetic moments are integer. 
The corresponding total spin magnetic moments for LSD, 
scenario SIC(l) and scenario SIC(2) are equal to 4.0 /is 
per formula unit, while for the third scenario the total 
spin moment is 2.0 /is per formula unit. The latter value 
results from the fact that for divalent ions, one additional 
d electron becomes localized. The type-resolved octahe- 
dral sites spin magnetic moments are substantially larger 
for the trivalent ions than for the divalent ions in all the 
SIC-LSD calculations. The spin magnetic moments of 
the tetrahedral Fe's are of similar magnitude but oppo- 
site sign to the trivalent octahedral ions. As there are 
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TABLE I: Total and type-decomposed spin magnetic mo- 
ments (in Bohr magnetons /is) for LSD and three differ- 
ent SIC-LSD scenarios for the cubic structure. Concern- 
ing types, only different Fe types are listed. Note that in 
the LSD, SIC(2) and SIC(3) calculations there is only one 
type of octahedral sites (B1=B2=B). Also given are the to- 
tal energy differences (in mRy/formula unit) with respect 
to the ground state configuration. These energy differences 
are further decomposed into the energy loss in band forma- 
tion/hybridization (AE_b_f) and the energy gain due to local- 
ization (AEs/c), both relatively to the respective energies of 
the ground state configuration, which is the all Fe''"'" scenario. 



Scenario Mtotai 








AE AEi3F 


AEsic 


LSD 4.00 


3.40 


-3.57 


-3.57 


894 -118 


1012 


SIC (1) 4.00 


4.00 


-3.57 


-4.08 


113 168 


-55 


SIC (2) 4.00 


4.02 


-3.90 


-3.90 








SIC (3) 2.00 


4.01 


-3.47 


-3.47 


374 475 


-101 



twice as many octahedral sites as tetrahedral ones, the 
total spin magnetic moment of magnetite in all the differ- 
ent scenarios is mostly due to the octahedral sites. The 
SIC-LSD spin magnetic moments of the oxygens are of 
the order of 0.02-0.30 ^b, depending on the scenario, 
while the spin magnetic moments on the empty spheres 
are still two orders of magnitude smaller. It is interesting 
that for the scenario with all octahedral ions in divalent 
configuration all oxygen moments are oppositely polar- 
ized to the octahedral Fe moments and are of the order of 
0.2-0.3 fiB- This of course is reflected in the much smaller 
total spin moment of this scenario. For the first and sec- 
ond scenarios the O moments are smaller, between 0.04 
and 0.17 /is, and polarized parallel to the octahedral Fe 
ions. The oxygen moments calculated within LSD are 
of the order of 0.1 /is and alligned anti-parallel to the 
octahedral Fe spin magnetic moments. 

Turning to Table II (orthorhombic structure), one can 
see similar behaviour to that of Table I (cubic struc- 
ture). Except for the LSD, the respective energy differ- 
ences for different scenarios are slightly smaller, but oth- 
erwise comparable to those for the cubic structure. Re- 
garding the type of solution, metallic or otherwise, again 
the scenario SIC(l) is insulating, for both orthorhom- 
bic structures, but the second (all trivalent) scenario for 
Pmc2i symmetry and LSD for the P2/c symmetry are 
not 100% spin polarized, as reflected by their respective 
total spin magnetic moments. Most importantly, how- 
ever, for the refined atom positions, as reflected by the 
structure with P2/c symmetry,^ the all trivalent scenario 
is half-metallic. Also, for the LSD solution in this struc- 
ture we find only a minute deviation from 100% spin po- 
larization with the DOS at the Fermi level of about 0.24 
states per formula unit and the total spin magnetic mo- 
ment of almost 4.00 fiB (with the accuracy better than 
0.5%). Compared to the other orthorhombic structure^ 
it appears that the refined atom positions are crucial for 
the half-metallic solution in the second scenario. Re- 
garding the type-resolved spin magnetic moments, again 



TABLE II: Total and type-decomposed spin magnetic mo- 
ments (in Bohr magnetons ^b) for LSD and three different 
SIC-LSD scenarios for the orthorhombic structures. Concern- 
ing types, only different Fe types are listed. In the LSD col- 
umn for the P2/c symmetry, 4.00(-) means that the magnetic 
moment is a tiny fraction less, but within better than 0.5% 
equal to 4.00 hb- Also given are the total energy differences 
(in mRy/formula unit) with respect to the ground state con- 
figuration. 





LSD 


SIC (1) 


SIC (2) SIC (3) 


Pmc2i symmetry 




4.00 


4.00 


4.02 


2.00 




3.39 


4.13 


4.14 


4.14 




3.32 


4.13 


4.13 


4.15 




3.30 


4.14 


4.16 


4.17 




3.34 


4.15 


4.16 


4.16 




-3.40 


-3.58 


-3.88 


-3.50 




-3.46 


-4.14 


-3.93 


-3.51 




-3.39 


-4.11 


-3.90 


-3.54 


Mf.b4 


-3.43 


-3.58 


-3.98 


-3.43 


Mf.b5 


-3.41 


-4.13 


-3.95 


-3.46 


MFese 


-3.37 


-4.11 


-3.94 


-3.48 


AE 


966 


105 





327 


P2/c symmetry 




4.00(-) 


4.00 


4.00 


2.00 




3.48 


4.08 


4.09 


4.08 




3.46 


4.06 


4.08 


4.07 


MFesi 


-3.53 


-3.57 


3.83 


-3.50 


MFes2 


-3.51 


-3.56 


3.83 


-3.51 




-3.67 


-3.56 


3.97 


-3.45 


MFeB4 


-3.39 


-3.55 


3.89 


-3.46 


MFeB5 


-3.40 


-4.01 


3.84 


-3.47 


MFese 


-3.57 


-4.05 


3.87 


-3.50 


AE 


841 


108 





337 



TABLE III: Listed are the numbers of spin-up and down va- 
lence electrons on the respective Fe-atom types for the LSD 
(two types only) and three different SIC-LSD calculations for 
the cubic structure. 

Scenario Fca^ Fea^ Fbbi^ Febi^ FeB2^ FeB2^ 

LSD 5.06 1.67 2.15 5.71 2.15 5.71 

SIC (1) 5.34 1.34 2.24 5.81 1.82 5.90 

SIC (2) 5.35 1.33 1.97 5.87 1.97 5.87 

SIC (3) 5.36 1.35 2.26 5.73 2.26 5.73 



we see that the trivalent tetrahedral and octahedral sites 
give rise to comparable values, with substantially smaller 
moments of the divalent Fe octahedral sites. Also, the 
oxygen spin moments are similar to those of the cubic 
structure. 

The spin resolved charges for various Fe atoms for all 
scenarios and structures are given in Tables III and IV. 
The first thing to note is that except for the LSD scenario 
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TABLE IV: Listed are the numbers of spin-up and down va- 
lence electrons on the respective Fe-atom types for the LSD 
(two types only) and three different SIC-LSD calculations for 
the orthorhombic structures. 





LSD 


SIC (1) 


SIC (2) SIC (3) 


Pmc2i symmetry 




5.41 


5.73 


5.74 


5.75 


FeAi^ 


2.02 


1.60 


1.60 


1.60 


FeA2^ 


5.35 


5.71 


5.71 


5.73 


FeA2^ 


2.03 


1.58 


1.58 


1.58 


FeA3^ 


5.33 


5.71 


5.72 


5.74 


FeA3^ 


2.03 


1.57 


1.56 


1.57 


FCAi^ 


5.28 


5.64 


5.65 


5.66 


Feaa^ 


1.94 


1.49 


1.49 


1.50 


Fesi^ 


2.28 


2.29 


2.04 


2.29 


FCBI^ 


5.68 


5.87 


5.92 


5.79 


FeB2^ 


2.21 


1.80 


1.97 


2.26 


FeB2^ 


5.67 


5.94 


5.90 


5.77 


FcB-i^ 


2.19 


1.77 


1.95 


2.20 


Fbb-j.^ 


5.58 


5.88 


5.85 


5.74 


FcBi^ 


2.20 


2.24 


1.92 


2.28 


FCBi^ 


5.63 


5.82 


5.90 


5.71 


FeB5^ 


2.19 


1.76 


1.91 


2.24 


FeB5^ 


5.60 


5.89 


5.86 


5.70 


FcBd^ 


2.22 


1.79 


1.93 


2.24 


FeBfi^ 


5.59 


5.90 


5.87 


5.72 


P2/c symmetry 


Feai^ 


5.05 


5.32 


5.33 


5.34 


Feai^ 


1.57 


1.24 


1.24 


1.26 


FeA2^ 


5.04 


5.31 


5.32 


5.33 


FeA2^ 


1.58 


1.25 


1.24 


1.26 


FCBI^ 


2.01 


2.08 


1.86 


2.09 


FCBl^ 


5.54 


5.65 


5.69 


5.59 


FeB2'^ 


2.04 


2.09 


1.88 


2.10 


FeB2^ 


5.55 


5.65 


5.71 


5.61 


FeB3^ 


1.88 


2.04 


1.72 


2.07 


FcB-i^ 


5.55 


5.60 


5.69 


5.52 


Feba^ 


1.92 


1.95 


1.67 


1.97 


Feba^ 


5.31 


5.50 


5.56 


5.43 


FCBS^ 


1.93 


1.57 


1.72 


1.97 


FeB5^ 


5.33 


5.58 


5.56 


5.44 


FeBG^ 


1.95 


1.65 


1.79 


2.05 


FeB6^ 


5.52 


5.70 


5.66 


5.55 



the electron charges of both spin components of the tetra- 
hedral Fe's do not differ too much among various scenar- 
ios. The situation is different for the octahedral sites. 
Especially, in the Verwey charge ordered scenario the dif- 
ferent octahedral ions show similar numbers in one spin 
component, but differ by about 0.42 electron in the other 
spin component of the cubic structure, which is due to 
localization of the sixth d electron. So, while in the ionic 



picture the valency of the two octahedral sites (Bl and 
B2) differs by one, in terms of the charge disproportion- 
ation they are only 0.32 electron different. The situation 
is similar for the orthorhombic structures. For the struc- 
ture with the Pmc2i symmetry the Bl and B4 sites are 
divalent, while B2, B3, B5 and B6 sites are trivalent. In 
the other orthorhombic P2/c structure, the sites Bl to B4 
are divalent, while B5 and B6 sites are trivalent. Again, 
the charge of spin down components of all the octahedral 
sites is of the order of 5.8-5.9 electron, the spin up charges 
of the divalent sites are about 2.2-2.3 electron, while for 
the trivalent sites they are of the order of 1.8 electron 
for the Pmc2i symmetry. For the P2/c symmetry the 
spin down charges of the octahedral sites fall between 
5.5 to 5.7 electron, while the spin down charges are of 
the order of 2.0 and 1.6-1.7 electron for the divalent and 
trivalent octahedral ions respectively. Consequently, the 
total charge disproportionation between the divalent and 
trivalent sites is up to 0.5 and 0.6 electron, respectively 
for the Pmc2i and P2/c symmetries. There is a spread 
in charges over the octahedral sites of the orthorhombic 
structures due to the reduced symmetry of these Fe sites. 
In the orthorhombic structure we also see charge dispro- 
portionation in the other three scenarios, namely LSD, 
scenario SIC(2) and SIC(3) - this of course does not oc- 
cur in the cubic phase. It is important to note that for 
the ground state configuration of all trivalent Fe's we find 
charge disproportionations of 0.2 and 0.36 electron for the 
Pmc2i and P2/c symmetries respectively. These calcu- 
lated charge disproportionations are in better agreement 
with the recent experimental study by Wright et al J than 
the ones of the Verwey charge ordered state. 



V. CONCLUSIONS 

In summary, in the present application of the first- 
principles SIC-LSD method to magnetite we have stud- 
ied three different arrangements of Fe's on the octahedral 
sites, namely the simple Verwey order, the case where, 
like the tetrahedral sites, all the octahedral sites have 
been occupied by Fe^^, and finally the scenario with all 
the octahedral sites occupied by the Fe^+ ions. For all 
these scenarios the tetrahedral sites have been occupied 
only by Fe'^^. Our total energy calculations, both for 
the cubic and orthorhombic structures, have shown that 
the second scenario, with all interstitial sites occupied by 
Fe'^"'", is the ground state solution, followed by the case 
of the Verwey charge order. In the latter case, we have 
calculated magnetite to be insulating with a gap of about 
0.35 eV ('^O.l eV in the orthorhombic phase), while the 
Fe'^"'" scenario, depending on structure, has given rise to 
half-metallic or metallic (though still with high spin po- 
larization) state. The insulating state for the Verwey sce- 
nario indicates that in the SIC-LSD description the O 2p 
bands are filled, as would be the case for the ionic model. 
However, no ionic charge distribution is obtained as sub- 
stantial Fe and O bonding and hybridization still occur. 
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Localizing one less Fe d electron (second scenario) results 
in disappearance of the insulating state and emergence of 
the half-metallic state, which is a non-trivial result, not 
predicted by the ionic model. 

The ground state scenario of all trivalent Fe's gives a 
charge disproportionation which is approximately half of 
what is found in the Verwey charge ordered state in the 
orthorhombic structure. The values of the ground state 
scenario are in line with the recent experimentsiiLiS 
The charge disproportionation in the ground state is sim- 
ilar to the one obtained from the LSD which would in- 
dicate that these charge disproportionations are not of 
an electronic origin but are determined by the structure. 
Namely, the inequivalent Fe positions give rise to differ- 
ent charge distributions, whilst the nominal Fe valency 
remains unaltered. 

Obviously we have managed to study only the three 
simplest charge ordering scenarios, and within these 
restrictions our results indicate the energetically un- 
favourable character of Fe^"*". We also find the same 
ground state configuration of all trivalent Fe in both the 
cubic structure and the distorted orthorhombic struc- 
tures. The calculations indicate no minority spin electron 
localization in the distorted structure, only a marginal 
decrease, in this structure, with respect to the cubic one, 
of 5% or so, in the energy difference between the ground 
state and the charge ordered one. This implies that a 
simple charge ordered state as described by scenario 1 
(SIC(l)) is not realized below the Verwey temperature. 

Our study has shown that the SIC-LSD is capable to 
study any kind of static charge order, although realiza- 
tion of more complicated charge arrangements, which 



might be of relevance to future studies, may be limited 
by memory and CPU of present computers. Given these 
limitations that we only studied a simple Verwey charge 
ordered state, we have concluded that the charge dispro- 
portionation below the Verwey temperature is structural 
in origin and all Fe ions occur in a trivalent state. This 
does not exclude the possibility that Fe304 below the 
Verwey temperature is described by a much more com- 
plicated charge ordered state. 
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FIG. 1: Spin-polarized total (top row), Fe-(middle row) and 0-(bottom row) partial densities of states from the LSD calculation 
(left column) in comparison with the SIC-LSD (right column) counterparts for the Verwey charge ordered phase in the inverse 
spinel structure. The partial DOS for empty spheres are of no significance and therefore have not been shown. 
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FIG. 2: Spin-polarized densities of states for the A-(top row) and Bl(middle row) -Fe sites in the LSD calculation (left column) 
(note that in the LSD there is only one type of octahedral sites (B1=B2=B)) in comparison to the respective A-(top row), 
Bl-(middle row) and B2-(bottom row) Fe sites from the SIC-LSD calculation (right column) in the Verwey charge ordered 
phase. 
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FIG. 3: Spin-polarized total (top), Fe- (middle), and O- 
partial (bottom) densities of states for the second charge order 
scenario, namely with trivalent configuration on all Fe-sites in 
the SIC-LSD implementation. 
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FIG. 4: Spin-polarized densities of states for the A- (top) and 
Bl- (bottom) Fe sites for the second scenario, namely with 
trivalent configuration on all, tetrahedral and octahedral, Fe- 
sites in the SIC-LSD implementation. Note that in this sce- 
nario there is only one type of octahedral sites (B1=B2=B). 



